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Numerical simulations of surface convection in a late M-dwarf 
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Abstract. Based on detailed 2D and 3D numerical radiation-hydrodynamics (RHD) simulations of time-dependent 
compressible convection, we have studied the dynamics and thermal structure of the convective surface layers of 
a prototypical late-type M-dwarf (T e s ~ 2800 K, log<? = 5.0, solar chemical composition). The RHD models 
predict stellar granulation qualitatively similar to the familiar solar pattern. Quantitatively, the granular cells 
show a convective turn-over time scale of ~ 100 s, and a horizontal scale of 80 km; the relative intensity contrast 
of the granular pattern amounts to 1.1%, and root-mean-square vertical velocities reach 240 m/s at maximum. 
Deviations from radiative equilibrium in the higher, formally convectively stable atmospheric layers are found to 
be insignificant allowing a reliable modeling of the atmosphere with ID standard model atmospheres. A mixing- 
length parameter of qmlt=2.1 provides the best representation of the average thermal structure of the RHD 
model atmosphere while alternative values are found when fitting the asymptotic entropy encountered in deeper 
layers of the stellar envelope (o:mlt = 1.5), or when matching the vertical velocity field (q:mlt=3.5). The close cor- 
respondence between RHD and standard model atmospheres implies that presently existing discrepancies between 
observed and predicted stellar colors in the M-dwarf regime cannot be traced back to an inadequate treatment of 
convection in the ID standard models. The RHD models predict a modest extension of the convectively mixed 
region beyond the formal Schwarzschild stability boundary which provides hints for the distribution of dust grains 
in cooler (brown dwarf) atmospheres. 
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1. Introduction 

Late-type M-dwarfs are fully convective stars where the 
convective flows penetrate far into the atmospheres reach 
ing o p tical depths as low a s 10~ 3 ( Allard fc Hauschildt 



1995). Allard et al. (1997) have reviewed the physical, 



spectroscopic, and photometric properties of these ob- 
jects. In the past, model atmospheres have typically failed 
to reproduce their spectroscopic and photometric 



proper- 



ties in I two respects: i) the ncar-IR spectral distribution 



chromospheric heating and/or the treatment of convection 
based on mixing- length theory (MLT). 

Hydrodynamical simulations of solar and stellar gran- 
ulation including a realistic description of radiative trans- 
fer have become an increasingly powerful and handy in- 
strument for studying the influence of convective flows 
on the the structure of late-type stellar atmospheres as 
well as on the formation of their spectra (e.g. Nordlund & 
Dravins 19901; ISteffen & Freytag 19911; iLudwig et al. 19941; 



(JHK. colors) where ; independent of the source of water 



m A^ a used, ™H*.i a ,11 agreed predict aa uader= 2000 D - Hith erto, model calculations have been exclusively 



vapor 

luminous if-band (relative to J), and ii) the optical My 
vs V — I color-magnitude relation, where all models sys- 
tematically predict bluer colors (i.e. being overluminous 
in V) than observed. 

Brett (1995) raised the possibility that this near-IR 
problem was due to models being "too cool in the upper 
photospheric layers", and suggested two possible causes: 



Freytag et al. 1996; Stein & Nordlund 1998; Asplund et al 
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performed for atmospheres where atomic lines are dom- 
inating the line blanketing. A possible next step in the 
development of hydrodynamical models is towards cooler 
atmospheres where molecular absorption dominates the 
atmospheric energy balance. Constructing hydrodynami- 
cal model atmospheres for cooler stars can shed light on 
the presently pressing shortcomings of the classical mod- 
els mentioned above. Regarding the considerable improve- 
ments in the quality of the molecular opacities and related 
atmospheric models, it becomes more and more important 
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to determine whether the treatment of convection by MLT 
is at the origin of the observed discrepancies. 

The basic questions we want to answer in this theoret- 
ical investigation are: Is mixing-length theory adequate 
to handle convection in the atmospheres of M-dwarfs? 
And if so, which mixing-length parameter «mlt is nec- 
essary to reproduce the various thermal and dynamical 
properties of an atmosphere (temperature profile in the 
line forming region, surface boundary condition connect- 
ing to stellar evolution models, convective velocities)? We 
start by describing some methodological aspects and the 
applied computer codes, in particular discuss the critical 
question of how accurately we can describe the complex 
radiative transfer within the hydrodynamical simulations. 
We continue by presenting our results which give some in- 
sight in what granulation looks like on the surface of an 
M-dwarf. We proceed with quantitative estimates of the 
mixing-length parameter, and discuss the consequences for 
conventional atmosphere modeling. Finally, we extrapo- 
late slightly beyond the existing hydrodynamical models 
proper, and suggest a scenario for the transport of dust 
grains in brown dwarf atmospheres due to convective over- 
shoot which is motivated from our present simulations at 
hotter temperatures. Often we refer to the Sun as our 
benchmark for comparison and assume some familiarity 
with its atmospheric properties. 

2. Methodological aspects, computer codes 

The aim is to model the atmospheric structure of a pro- 
totypical late M-dwarfs as realistically as possible, with 
a focus on the interplay of convective flows and radiative 
transfer. Being well aware of the limitations in our models, 
we took, whenever possible, a differential approach in try- 
ing to reduce the influence of systematic uncertainties on 
the outcome of our investigation. This concerned mostly 
the dimensionality of the problem: the multi-dimensional, 
time dependent approach adopted in the hydrodynamical 
simulations versus the one-dimensional, time independent 
approach adopted in classical stellar atmospheres. We en- 
sured that the numerical treatment (i.e. implemented mi- 
crophysics, representation of radiative energy transport) 
in the two "worlds" was as similar as possible. We em- 
ployed various computer codes whose names and main 
characteristics we introduce below. We elaborate on spe- 
cific aspects critical for the investigation in more detail 
later. 

RHD: A radiation- hydrodynamics code developed by 
A. Nordlund and R.F. Stein (see (Btein fc Nordlund 1998| , 
and references therein) for modeling stellar atmospheres 
in two or three spatial dimensions. It implements a con- 
sistent treatment of compressible gas flows together with 
non-local radiative energy exchange. The radiative trans- 
fer is treated in LTE approximation, the wavelength de- 
pendence of the radiation field is represented by a small 
number of wavelength bins (see below). Open lower and 
upper boundaries, as well as periodic lateral boundaries 
are assumed. The effective temperature of a model (i.e. the 



average emergent radiative flux) is controlled indirectly by 
prescribing the entropy of inflowing material at the lower 
boundary. Magnetic fields are neglected. 

LHD: A ID Lagrangian hydrodynamics code devel- 
oped by one of the authors (HGL) used to calculate stan- 
dard stellar atmospheres which can be compared with 
results obtained with RHD. Besides the reduced dimen- 
sionality, the adopted physics (equation of state, radiative 
transfer scheme) is the same as in RHD. The convective 
energy transport is based on MLT. In this paper we adopt 
the MLT formulation of |Mihalas (1978| ). Excluding one 
exception, all values of the mixing-length parameter are 
given with reference to Mihalas' formulation. 

PHOENIX: A ID model atmosphere code developed 
by two of the authors (PHH & FA, for a detailed descrip- 
tion see Hauschildt & Baron 1999). It implements a treat- 
ment of the wavelength dependence of the radiation field 
with high resolution based on direct opacity sampling. In 
this investigation PHOENIX served as opacity data base, 
and was used for assessing the quality of the simplified ra- 
diative transfer employed in the hydrodynamics codes. 

ATLAS6: A version of the ID model atmosphere code 
developed by R. Kurucz QKurucz 1979, ). It served as addi- 
tional opacity data base. 

Table [J summarizes the properties of the hydrodynam- 
ical models discussed in the paper. Model A-3D is our 
M-dwarf reference model. The twice as large model B-3D 
was primarily calculated for checking effects of the do- 
main size, the solar model S-3D was added for assessing 
scaling properties with effective temperature and gravita- 
tional acceleration. We note that the solar model is not 
strictly differentially comparable to the M-dwarf models 
since it is employing different opacity sources and equation 
of state which stem from the U ppsala stellar atmosphere 
package ( Gustafsson et al. 1975 ). We do not consider this 
as particularly critical since the physical conditions in the 
atmospheres of M-dwarfs and the Sun are so different that 
one looses the advantages of a differential approach any- 
way. The 2D models C-2D and D-2D are models which 
were considered in the forefield to investigate effects of the 
numerical resolution. They are based on ATLAS6 opac- 
ities without contributions of molecular lines and employ 
grey radiative transfer. Due to the different input physics 
their behavior is qualitatively different from the more re- 
alistic 3D models. Despite their shortcomings, they are 
of interest for qualitatively understanding the interplay 
of convection and radiative transfer in optically thin re- 
gions, and therefore will be discussed in more detail. As 
is apparent from the table, all hydrodynamical M-dwarf 
models show very small fluctuations around their aver- 
age effective temperature. This reflects the fact that the 
horizontal and temporal fluctuations of all quantities are 
small compared to the Sun — a central feature of the at- 
mospheres of late- type M-dwarfs. 
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Table 1. The hydrodynamical models discussed in the paper: Name is used to identify a model, Dim. the dimen- 
sionality of the model, Mesh the number of grid points (XxYxZ, where Z denotes the vertical, X and Y the horizontal 
directions), Size the geometrical size of the computational domain, Opacities the employed opacity source, -ZVobm 
the number of OBM bins ("Opacity Binning Method", for a description see Sect. 2.2) for calculating the radiative 
transfer, T c g the resulting effective temperature of the model including an estimate of its RMS fluctuations, \ogg the 
gravitational acceleration, and Modelcode an internal model identifier. 



Name 


Dim. 


Mesh 


Size [km] 


Opacities 


Aobm 


Teff [K] 


logs 


Modelcode 


A-3D 


3 


125x125x82 


250x250x87 


PHOENIX 


4 


2789.1 ±0.7 


5.0 


d3gt30g50nl8 


B-3D 


3 


250x250x82 


500x500x116 


PHOENIX 


4 


2799.8 ± 0.2 


5.0 


d3gt30g50nl9 


S-3D 


3 


125x125x82 


6.0x6.0x3.2 (Mm) 


Uppsala 


4 


5640 ± 14 


4.44 


sun3d 


C-2D 


2 


125x82 


250x80.3 


ATLAS6 


1 


2774 ± 1.6 


5.0 


d2gt30g50n8 


D-2D 


2 


251x162 


250x80.3 


ATLAS6 


1 


2770 ± 1.6 


5.0 


d2gt30g50n9 




4 5 6 

log 10 P gas [dyn/cm z ] 

Fig. 1. Comparison of the specific heat at constant pres- 
sure from various equations of state in a representative 
M-dwarf: RHD (solid), and PHOENIX (dash-dotted). 
If one artificially removes the contribution of H2 molecu- 
lar formation c p is largely reduced (dotted). For further 
comparison the run of c p in a solar model atmosphere 
is depicted (dashed). For the M-dwarf model a logarith- 
mic Rosseland optical depth scale is indicated by the tick 
marks near the abscissa. 



2.1. The equation of state 

Figure |] shows that the equation of state (EOS) which 
is employed in the hydrodynamical simulations and LHD 
models is very similar to the PHOENIX EOS. We do not 
expect significant systematic differences by applying these 
two different equations of state in our various model cal- 
culations. The equations of state of RHD and PHOENIX 
treat ionization and molecular formation assuming Saha- 
Boltzmann statistics. Since non-ideal effects arc not pro- 
nounced in atmospheres of M-dwarfs the inclusion of H2 
molecular formation is the main ingredient required to ob- 
tain a realistic description of the thermodynamics of the 
stellar plasma. 

Indeed, Fig. |l| demonstrates that in the M-dwarf at- 
mosphere H 2 molecular formation is the most important 
contributor for increasing the specific heat above the value 



associated with the purely translatorial degrees of free- 
dom. In the Sun the dominant contributor is the hydrogen 
ionization which has an even more dramatic effect on the 
heat capacity of the stellar plasma. I.e. from the perspec- 
tive of the content of latent heat in the gas flows M-dwarf 
atmospheres are not particularly extreme objects. 

2.2. Radiative transfer 

An important problem when modeling M-dwarf envelopes 
is the treatment of the large number of absorption lines 
in their atmospheres which are mostly of molecular origin. 
The complex wavelength dependence of the radiation field 
is illustrated in Fig. [|. While it is already a formidable task 
to treat the frequency dependence of the radiation field 
in ID model atmospheres, this is even more the case in 
hydrodynamical models where one has to account for the 
3D geometry of the flow field and its temporal evolution. 
Present computer capacity allows only for a very modest 
number of frequency points to be included the modeling 
of the radiation field within a hydrodynamical simulation. 
However, the situation is somewhat alleviated. For the 
interaction of hydrodynamics and radiative transfer only 
the frequency integrated net amount of radiative heating 
(or cooling if negative) 



Qrad = 47rp / Xv{Jv ~ S v ) dv 



(1) 



(p denotes the density, \v the monochromatic absorp- 
tion coefficient, J„ the angular average of the intensity, 
S v the source function) is of relevance. I.e. we do not 
need to know the detailed frequency dependence of the 
radiative heating. Perhaps surprisingly, one can obtain 
a reasonably accurate description of the radiative heat- 
ing with very few frequency points by a judicious choice 
of the frequency discretization: frequencies with a simi- 
lar run of monochromatic optical depth are grouped to- 
gether into opacity "bins" . In the context of hydrody- 
namical stellar atmospheres this Op acity Binning M ethod 
( OBM) was originally introduced by Nordlund (1982 ), and 
some variants and refinements were subsequently devel- 
oped (Ludwig et al. 1994). The basic idea is to group fre- 
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Fig. 2. Scatter plot of standard optical depth where the 
monochromatic optical depth reaches unity as a func- 
tion of wavelength in a PHOENIX model atmosphere 
with T ff =2900 K, logg=5.3, and solar chemical compo- 
sition. Standard optical depth is the optical depth scale 
in the continuum plus opacity contributions of H2 and 
H2O molecular lines at A = 1.2 /im. The horizontal lines 
mark the thresholds which were used for defining the OBM 
opacity bins (see text). The curve in the lower part of the 
panel is the emergent flux (arbitrarily scaled) of the model. 



quencies which reach monochromatic optical depth unity 
within a certain depth range of the atmosphere into one 
frequency bin. Effectively, this combines frequencies into 
the same bin which belong to the continuum, weaker, and 
stronger spectral lines. The opacity in each bin is an aver- 
age (Rosseland or Planck) over the frequency range which 
is represented by that bin. With these average opacities, 
supplemented by the frequency integrated source function, 
the equation of radiative transfer has to be solved only 
once per bin to obtain the formal solution of the radiative 
transfer for a given flow field. The binning procedure is il- 
lustrated in Fig. H where the horizontal lines separate the 
depth ranges — as measured by a frequency independent 
standard optical depth — which are sorted into the vari- 
ous OBM bins. In the present investigation we used 4 bins 
in our non-grey hydrodynamical models, representing the 
continuum and line wings, as well as the cores of weak, 
medium strong, and strong spectral lines. 

It is clear that the sorting procedure described before is 
specific to the stellar atmosphere under consideration and 
has to be repeated as soon as the atmospheric parameters 
differ widely from the one where the sorting was done. 
We used a model atmosphere calculated with PHOENIX 
at T e ff = 2900 K and log<? = 5.3 as reference atmosphere 
for the grouping. It is sufficiently close to the atmospheric 
parameters of the hydrodynamical models. This has been 
checked by studying the performance of the sorting when 
applied to differing atmospheric parameters (log g down 
to 3.0, and T cff up to 3300 K). 

The sorting criterion that t\ = 1 should fall within 
a certain depth range in the atmosphere does not guar- 



antee that the overall depth dependence of t\ is similar 
for all frequencies grouped together — a condition for 
allowing the interchange of the solution of the transfer 
equation with the frequency integration when evaluating 
Q ra d. In particular, the simultaneous presence of atomic 
and molecular lines can lead to significantly different func- 
tional forms of the monochromatic optical depth at differ- 
ent frequencies: deeper regions of the atmosphere too hot 
to allow for molecule formation might be dominated by 
atomic lines while higher and cooler regions which allow 
for molecule formation might be dominated by molecular 
lines. If the atomic and molecular lines emerge from differ- 
ent elements there is no physical connection between them 
leading to an uncorrelated behavior in optical depth of 
deeper and higher layers with frequency. Such a situation 
would be unfavorable for the OBM. However, the OBM 
is rather successful in reproducing the heat exchange be- 
tween radiation and matter — as evident from Fig. p[ This 
is linked to the statistical dominance of molecular absorp- 
tion in all radiative layers of the rather cool atmosphere 
under consideration. 

Another point concerning the present implementation 
of the OBM is our usage of global Rosseland means — i.e. 
Rosseland averages over the whole frequency range — for 
representing the average opacity in each bin. For the the 
continuum we took the Rosseland means themselves while 
we scaled them by factors of 10 1 , 10 2 , and 10 3 for the 
bins representing the successively stronger lines. The scal- 
ing factors correspond to the factors of 10 in optical depth 
selected as thresholds for the sorting procedure which are 
log 10 T s tandard = -0.5,-1.5, and - 2.5 (see Fig. ||). The 
basic assumption behind this approach is that the line 
opacity scales with temperature and pressure like the con- 
tinuous opacity. The source function has been integrated 
over the frequency ranges of the individual bins. While it 
is certainly not the optimal representation of the opaci- 
ties it was dictated by the lack of a detailed tabulation 
of the monochromatic opacities over the full temperature- 
pressure range of interest. Work is presently under way to 
generate such tabulation which is a non-trivial task due to 
the enormous amount of line data which need to be pro- 
cessed. Considering the various approximations described 
before one might ask whether the OBM is a real improve- 
ment beyond a grey description. Eventually, this can be 
tested quantitatively by comparing ID atmospheres com- 
puted with high frequency resolution or the OBM. Such 
a comparison serves as ultimate indicator of the perfor- 
mance of the OBM. 

Figure ^| shows a comparison of ID model atmospheres 
calculated with PHOENIX employing opacity sampling 
with a wavelength resolution of » 2 A, and LHD employ- 
ing the OBM as described before. Comparing flux con- 
stant models in radiative and radiative-convective equilib- 
rium shows a good correspondence of the resulting equi- 
librium temperature profiles. While the radiative equilib- 
rium models are less important for the investigation of M- 
dwarf atmospheres, they were added to the comparison 
to show the similarity in the radiative transport proper- 



Ludwig et al.: Numerical simulations of surface convection in a late M-dwarf 




5 6 
log 10 P gas [dyn/cnf] 

Fig. 3. Comparison of ID standard model atmospheres 
(T e ff=2800K, log 5=5.0) in purely radiative (dashed 
lines) and radiative-convective (q;mlt = 1.0, solid and 
dotted lines) equilibrium. The thick lines are based on 
LHD atmospheres using the approximate OBM radiative 
transfer, the thin lines are based on PHOENIX models 
using opacity sampling with high wavelength resolution. 
The dotted line depicts an LHD model employing grey 
radiative transfer. A logarithmic Rosseland optical depth 
scale is indicated by the tick marks near the abscissa. 



4 5 
log ,0 P gas [dyn/cm 2 ] 

Fig. 4. Various characteristic time scales as a function 
of pressure in an M-dwarf model: Brunt- Vaisala period 
(solid), Kelvin- Helmholtz time scale (dashed), radia- 
tive relaxation time based on ATLAS6 Rosseland opaci- 
ties (dotted), PHOENIX Rosseland (dashed-dotted), 
and PHOENIX Planck (triple-dot-dashed) opacities. 
A logarithmic Rosseland optical depth scale is indicated 
by the tick marks near the abscissa. The kink in the run 
of the Brunt- Vaisala period is located at the boundary 
between the convectively stable and unstable part of the 
model. 



ties independent of influences of the convective transport. 
More important are the models in radiative-convective 
equilibrium since they are closer to the actual physical 
situation. For judging the correspondence between opac- 
ity sampling and OBM models one should take the grey 
model as benchmark which represents a strongly simplified 
(equivalent to one frequency bin) description of radiative 
transfer. Clearly, the OBM matches much more closely the 
PHOENIX opacity sampling model. The cooling of the 
higher atmosphere by lines is reasonably represented as 
well as the backwarming of the deeper layers. 

The most important deviation between opacity sam- 
pling and OBM model occurs in the layers where the tran- 
sition from convectively to radiatively dominated energy 
transport takes place (around logP = 5.3). The OBM 
stratification becomes noticeably cooler. This was traced 
back to an insufficient heating of the gas in the OBM con- 
tinuum bin. This in turn is probably related to the use 
of Rosseland averages for the continuum opacity in these 
optically thin regions. Planck averages would be more suit- 
able. But again, at the time this work was performed only 
globally averaged opacities were available. Rosseland and 
Planck averages differ by at least a factor of 100 in these 
regions, making an ad-hoc switching from one to the other 
problematic. However, we consider the remaining devia- 
tions not as vital, in particular with respect to the cali- 
bration of the mixing- length parameter which we describe 
later in this paper. The calibration is a result of a differ- 
ential comparison of models which all base on the OBM. 

For completeness, we finally remark that only the 
exchange of energy was considered within the radiative 



transfer. The exchange of momentum was neglected; the 
prevailing relatively high mass densities combined with 
low radiative fluxes render radiation pressure unimportant 
for the structure of M-dwarf atmospheres. 

3. Characteristic time scales 

To get an overview of the problem Fig. ^ shows var- 
ious characteristic time scales in a representative M- 
dwarf model (LHD model with T c ff=2790K, log #=5.0, 
(*mlt=T-0). The expressions which we applied in the cal- 
culation of the time scales are summarized in the ap- 
pendix. The time scales have been evaluated under simpli- 
fying assumptions, and should therefore be taken as order 
of magnitude estimates only. 

The radiative time scales based on ATLAS6 and 
PHOENIX Rosseland mean opacities are rather similar in 
the deeper layers while being very different in the optically 
thin regime. This emphasizes the great influence of molec- 
ular absorption which becomes important at cooler tem- 
peratures and which is not included in the ATLAS6 opaci- 
ties. The radiative time scales from PHOENIX Rosseland 
and Planck opacities differ also significantly. Comparing 
the radiative time scales to a dynamical or convective time 
scale as given by the Brunt- Vaisala period leads us to ex- 
pect that one gets a qualitatively different behavior de- 
pending on the treatment of the radiative transfer. We 
will see that models based on a frequency-independent 
Rosseland opacities show significant deviations from ra- 
diative equilibrium conditions. A more realistic treatment 

in the optically thin regions more closely represented 
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by the Planck mean opacities — results in an atmospheric 
structure closer to radiative equilibrium. In all cases, one 
expects an almost adiabatic structure in the deeper layers 
since the radiative relaxation time is much longer than the 
convective time scale. 



3.1. A remark on thermal relaxation 

RHD uses an explicit numerical scheme for advancing 
the solution in time. It is well known that in an explicit 
scheme the time step is limited to a fraction of the dy- 
namical time scale, more precisely to the limit given by 
the Courant-Friedrichs-Levy stability condition. Together 
with the amount of available computer — and wallclock — 
time, this limits the time interval which can be covered by 
a 3D RHD model to about 10 4 s of stellar time. Figure f| 
seemingly implies that it is impossible to obtain a ther- 
mally relaxed hydrodynamical model within this time in- 
terval since the Kelvin- Hclmholtz time scale of the deeper 
layers is at least an order of magnitude larger. Contrary 
to the impression given in Fig. |], in the multi-dimensional 
hydrodynamical simulations this does not pose a problem 
since the thermal relaxation of the model is not governed 
by the time scale for the exchange of energy as expressed 
by the Kelvin-Helmholtz time. In fact, the thermal evolu- 
tion of the deeper layers of the RHD models is governed 
by the time scale for the exchange of mass in these layers 
which is much shorter (as quantified below, see Fig. |l6| ). 
Due to the exponential run in density of the atmosphere 
the mass exchange consists primarily of the replacement 
of mass by fresh material stemming from deeper layers. 
In the hydrodynamical model it is ultimately fed into the 
computational domain at the lower boundary. The con- 
vective energy flux is the net effect of the energy trans- 
ported by counteracting, opposing mass currents in and 
out of a test volume. The advective — as opposed to dif- 
fusive — nature of the mass transport leads to a situation 
where there can be a large imbalance in the energy con- 
tent of the mass currents: the in-coming mass elements 
carry the energy (strictly speaking the specific entropy) 
of much deeper layers (in the hydrodynamical model the 
entropy at the lower boundary), while the out-going el- 
ements just carry the local energy density. This implies 
an energy flux much larger than the one close to equi- 
librium conditions. The system is driven much faster to- 
wards equilibrium than implied by the classical Kelvin- 
Hclmholtz time scale, which assumes an energy flux as 
encountered close to equilibrium conditions. 

We have argued from the perspective of our hydrody- 
namical models. But the fast thermal relaxation is driven 
by a physical mechanism implying that a real stellar con- 
vection zone does not operate any differently. Hence, we 
want to stress that to our understanding the thermal re- 
laxation of a convective layer is generally not governed 
by the classical Kelvin-Hclmholtz time scale but by the 
usually shorter time scale of mass exchange. 



In the one-dimensional LHD models no mass exchange 
takes place and the convective energy flux is computed 
from MLT, in which convection is modeled as diffusion 
process. Under such circumstances the time scale of the 
thermal evolution of the convective layers is indeed com- 
parable to the Kelvin-Hclmholtz time scale as shown in 
Fig. ^. Of course, covering the time interval necessary for 
the slower thermal relaxation in a ID model run poses no 
problems due to the largely reduced computational costs. 

4. 2D precursor models 

Models C-2D and D-2D arc 2D models (see Tab. @) which 
were constructed to gain experience regarding the size 
of the computational domain and grid resolution. They 
are based on grey radiative transfer utilizing ATLAS6 
opacities which do not include contributions of molecular 
lines. The choice of the ATLAS 6 opacities was dictated 
by the lack of more realistic opacities during the initial 
phase of the project. From Fig. |] one might readily con- 
clude that their atmospheric structure will be dominated 
by convection since the radiative relaxation times in the 
atmosphere are long as compared to the dynamical time 
scale. Indeed, when starting from a temperature structure 
taken from a LHD model based on MLT, we find a rapid 
cooling of the originally radiatively stratified atmospheric 
layers by convective overshooting (see Fig. ||). Convection 
tends to transform the stratification into a purely adia- 
batic one since the radiative heating is too weak to keep 
the temperature close to the radiative equilibrium value. 
The transformation into an almost adiabatic stratification 
takes too long to be covered within the multi-dimensional 
hydrodynamical simulations. However, we conducted nu- 
merical experiments with LHD where an ad-hoc velocity 
field mimicking the convective overshoot was put into the 
atmospheric regions which are formally stable according 
to the Schwarzschild criterion. The LHD models indicate 
that the asymptotic temperature profile tends to a fully 
adiabatic stratification. 

Due to the unrealistic opacities, the models cannot 
give a good representation of a real M-dwarf atmosphere. 
They are nevertheless interesting from a numerical point 
of view. The resolution of the less resolved model C-2D 
is already sufficient to represent the convective transport 
properties, in particular in the important transition re- 
gion from convectively to radiatively dominated energy 
transport: after an initial relaxation phase we find in both 
models that the minimum temperature drops linearly with 
a rate of ^ = — 0.011 K/s. Looking at further diagnostics 
at comparable instances during the secular evolution of 
the model runs we observe that RMS velocities are simi- 
lar within a 20 % level. The higher resolved model shows 
somewhat more small-scale features, and its downflows 
appear more concentrated. However, we are primarily in- 
terested in the balance of convective to radiative energy 
transport. The temperature drop rate is a convenient mea- 
sure of the relative efficiency of both processes. Hence, we 
conclude from the similarity of the drop rates that the 
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Fig. 5. Temporal evolution of the average temperature 
profile of model C-2D in steps of 5ks showing a suc- 
cessive cooling of the layers around logP = 5.5 (solid 
lines). For comparison a LHD model with «mlt = 1-0 
is plotted (dash-dotted). All models are based on grey 
radiative transfer employing ATLAS6 opacities devoid of 
molecular contributions. All models have the same atmo- 
spheric parameters (T c ff=2770K, logg=5.0). A logarith- 
mic Rosseland optical depth scale is indicated by the tick 
marks near the abscissa. 

resolution of the less resolved model applied in 3D sim- 
ulations is sufficient to model the transport properties of 
the convective flows. 

Another conclusion which can be drawn from these 
models is that standard MLT models can give quite mis- 
leading predictions of the atmospheric temperature struc- 
ture if the radiative relaxation time is long in comparison 
to convective time scales. In other words, one should be 
cautious when the coupling of the temperature structure 
to the radiative equilibrium temperature is so weak that 
overshoot of low amplitude happens essentially adiabat- 
ically. In the following we shall see that M-dwarf atmo- 
spheres are unlikely to exhibit such conditions. 

5. 3D models: results and discussion 

5.1. General flow structure 

Figure ^ shows a typical snapshot of the emergent inten- 
sity during the temporal evolution of model B-3D. For 
comparison, Fig. |?] shows a similar snapshot from the so- 
lar run S-3D. Note, that the spatial and intensity scaling 
is very different in the reproductions, hence, only rela- 
tive geometrical properties should be directly compared. 
The average relative RMS intensity contrast of the gran- 
ular pattern amounts to 1.1% in the M-dwarf as opposed 
to 16% in the solar case. The first thing to note is that 
surface convection in an M-dwarf produces a granular pat- 
tern qualitatively resembling solar granulation: bright ex- 
tended regions of upwelling material which are surrounded 
by dark concentrated lanes of downflowing material. The 
dark lanes form an interconnected network. Looking more 
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Fig. 6. Typical snapshot of emergent intensity during the 
evolution of model B-3D. The intensity contrast amounts 
to <5/ rm s = 1.2 % at this instant in time. Note, that the con- 
trast of the intensity pattern in the graphics is enhanced 
with respect to its real appearance. 

closely, granules are less regularly delineated in M-dwarfs, 
the inter-granular lanes show a higher degree of variability 
in terms of their strength. A feature which is uncommon 
in the solar granulation pattern are the dark "knots" (e.g. 
at x = 100 km and y = 230 km in Fig. ^|) found in or 
attached to the inter-granular lanes. The knots are associ- 
ated with strong downdrafts which carry a significant ver- 
tical component of angular momentum. The width of the 
inter-granular lanes to the typical granular size is smaller 
in M-dwarfs. Inspecting the velocity field (not shown) in 
vicinity of the continuum forming layers shows less pro- 
nounced size differences. This indicates that the relatively 
broader lanes in the solar case are the result of a stronger 
smoothing of the temperature field due to a more intense 
radiative energy exchange, i.e. the effective Peclet num- 
ber of the flow is larger around optical depth unity in 
M-dwarfs. 

The different magnitude of the intensity contrast al- 
ready indicates that horizontal fluctuations of the ther- 
modynamic quantities are small in M-dwarf atmospheres. 
Figure |^ shows the run of the relative temperature and 
pressure fluctuations in model A-3D. Plotted are long term 
temporal and horizontal Q RMS averages. As we will argue 
later, the fluctuations in the higher atmosphere are likely 
to be overestimated in the model. But even taken at face 
value they are quite modest. The low level of fluctuations 

The horizontal averaging is performed on the geometrical 
depth scale. Due to the small horizontal fluctuations, alterna- 
tive choices (e.g. on the optical depth scale) give similar results. 
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Fig. 7. Like Fig. |, but for the solar model S-3D. The in- 
tensity contrast amounts to i5/ rms = 16.7% at this instant 
in time. Note, that for reproduction purposes spatial as 
well as intensity scale are altered in comparison to Fig. [| 



in the thermodynamic quantities is accompanied by small 
flow velocities, the maximum Mach number amounts to 
6.5 % in model A-3D. MLT models show that Mach num- 
bers drop to even lower values as one goes to lower effec- 
tive temperatures which are encountered in the regime of 
brown dwarfs. Our findings may have a direct bearing on 
the dust formation conditions in such objects. Our hydro- 
dynamical simulations support the view that dust form- 
ing layers in cool main sequence objects experience only 
small variations around their mean thermodynamic state. 
This view is clearly at odds with the scenario discussed 
by Helling ct al. (200l| ) who study the dust formation in 
turbulent brown dwarf atmospheres. Helling et al. assume 
thermodynamic fluctuations of order unity. 

Another important feature, which distinguishes M- 
dwarf atmospheric conditions from solar ones, is the ex- 
tent of the convective layers. As evident from Figures [ll] 
and [l2|, the convective motions reach much lower optical 
depth (log TRossciand ~ —1-5) in the M-dwarf than in the 
Sun (log TR OSSe iand ~ 0). Two factors contribute to this 
behavior. The temperature gradient in non-grey radiative 
equilibrium is significantly steeper in the M-dwarf than 
in the Sun, presumably due to efficient cooling of the at- 
mospheric layers by molecular lines. Furthermore, the H2 
molecule formation reduces the adiabatic gradient in the 
M-dwarf atmosphere, while in the Sun the hydrogen re- 
combination is essentially completed in subphotospheric 
layers. The different radiative and thermodynamic condi- 
tions favor the presence of convection in M-dwarf atmo- 
spheres. 
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Fig. 8. Relative horizontal temperature (solid) and pres- 
sure (dashed) fluctuations of model A-3D. The fluctua- 
tions in the region log P < 5 are likely to be overestimated 
in the model, i.e. should be considered as upper limits. 
The dotted line indicates the approximate location of the 
Schwarzschild boundary. A logarithmic Rosseland optical 
depth scale is indicated by the tick marks near the ab- 
scissa. 



In the deeper layers we observe the tendency — also 
known from solar simulations (see Stein & Nordlund 1998) 
— that the granular network of downflows decays into iso- 
lated downdrafts. Our models are rather shallow reaching 
only 2.3 pressure scale heights below optical depth unity. 
Thus we cannot follow the change of flow topology as far 
as has been done for solar models, but within the limited 
depth range comprised by our models we do not see in- 
dications of a qualitatively different behavior as found in 
the Sun. 



5.2. Horizontal scales 

In the following we shall discuss spatial power spectra of 
intensity and vertical velocity. Note, that in the Figures || 
and [l^ we display power per logarithmic wavenumber in- 
terval and not per unit wavenumber — the more common 
choice. This allows visually for a more direct identification 
of the power carrying scales. However (in humble respect 
of Kolmogorov's achievements) , we labeled the power laws 
drawn for comparison with the familiar spectral index of 
power per unit wavenumber. The spectra are temporal av- 
erages over one to a few convective turn-over times and 
many convective cells, so that they are statistically repre- 
sentative. 

Figure ^| shows a comparison of power spectra of the 
emergent intensity (more precisely: the intensity in the 
OBM continuum bin in vertical direction at the upper 
boundary of the computational volume) for the models A- 
3D, and B-3D as well as the solar model. Figure [l(] shows 
a corresponding comparison of the vertical velocity com- 
ponent measured at the layer where its RMS value reaches 
the maximum. To facilitate a comparison between the M- 
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the height of the maxima in the M-dwarf models, while the 
horizontal scale of the solar model had to be shrunk by a 
factor of 15. This factor closely corresponds to variation of 
the pressure scale height at the surface (factor 12.5) and 
follows the scaling found for hotter models by Freytag 
et al. (1997). Different from Freytag et al., we accounted 
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Fig. 9. Power spectrum of the emergent intensity pattern 
for model A-3D (solid), and B-3D (dashed), as well as 
the solar model S-3D (dotted). The spectrum of the so- 
lar model was scaled in power as well as wavenumber to 
match the peak of power in the M-dwarf models. Lines 
with slopes 4£ oc fc~ 5 / 3 and 4r oc fc -17 / 3 are shown for 
comparison with expectations from simplistic turbulence 
theory. The total relative intensity contrasts amount to 
5I lms = 1.12% (A-3D), 1.12% (B-3D), and 16.2% (S-3D). 
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Fig. 10. Power spectrum of the vertical velocity pattern, 
line styles as in Fig. |^. The solar power spectrum was 
scaled to match the peak of power of the M-dwarf mod- 
els. Thle velocity pattern was taken from the layer where 
the RMS vertical velocity reaches its maximum, with abso- 
lute values of w rmSi , nax = 236m/s (A-3D), 241m/s (B-3D), 
2600 m/s (S-3D). 



dwarf models and the solar model, the solar model was 
arbitrarily scaled in power and wavenumber so that the 
maxima and position of the power distributions matched. 
For the M-dwarf models the hump in power at the highest 
wavenumbers is likely an artifact of an imperfect choice 
of parameters controlling the small scale dissipation and 
should be ignored. 

Fluctuations in the Sun are an order of magnitude 
larger and had to be reduced accordingly in order to match 



for the change of the mean molecular weight of the stellar 
gas since the H2 molecular formation leads to a significant 
increase in M-dwarf atmospheres. The peak power of the 
M-dwarf models is located at an absolute scale of 80 km 
for intensity and 63 km for velocity structures. 

Comparing the width of the power distribution at 1 dex 
below peak level, the M-dwarf model B-3D displays a 
range of scales of 1.4 dex, the solar model of 1.1 dex. As 
is evident from the figures, the solar scale distribution is 
slightly but noticeably narrower. This might be traced 
back to the lower Peclet number of the plasma encoun- 
tered at the surface of an M-dwarfs as opposed to solar 
conditions. The formation of small scale structures is more 
strongly suppressed by the radiation field under solar con- 
ditions, leading to an overall narrower power distribution. 
Model A-3D apparently lacks some larger scale power indi- 
cating that the computational box is somewhat too small 
in this case. However, since the total fluctuations are very 
similar in the models we conclude that the relevant struc- 
tures are captured in model A-3D which serves as our 
standard model. 

Towards high wavenumbers the M-dwarf spectra do 
not show a clear power law behavior. In any case, they 
cannot be fitted with a —5/3 slope, and are at best roughly 
compatible with a —17/3 slope. We would like to empha- 
size that despite the very different physical parameters of 
M-dwarf atmospheres this is qualitatively what was al- 
ready found in solar models, and is at odds with expec- 
tations for homogeneous and isotropic turbulence. While 
smaller scales in the hydrodynamical models are definitely 
influenced by the artificial viscosity, it is still somewhat 
surprising that the asymptotic behavior in the Sun as well 
as in the M-dwarf is so similar despite their quite different 
characteristic flow speeds (i.e. Mach numbers). While our 
present models cannot make definite statements about the 
properties of small scale turbulence, a "non-Kohnogorov" 



behavior cannot be ruled out neither. See Nordlund et al 



1997 ) for an in-depth discussion of this issue in a solar 



context. 

5.3. Mean temperature stratification and 

corresponding effective mixing-length parameter 

Figures [ll], [l^, and [ll] show a comparison of the ther- 
mal structure of our "best" hydrodynamical model A- 
3D with standard ID atmosphere models computed with 
LHD. Again, we point out that the treatment of the mi- 
crophysics (equation of state and opacities) is equivalent 
between RHD and LHD models. Differences lie mainly 
in the dimensionality of the models — 3D versus ID - 
and the related MLT treatment of the convective energy 
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Fig. 11. Overview of the temperature-pressure relation of 
hydrodynamical model A-3D (solid) and three standard 
mixing-length model atmospheres (dash-dotted) with 
«mlt = 1-0, 1.5, 2.0. The «mlt = 1-0 model is the hottest 
at high pressure while it is coolest around logP = 5. A 
logarithmic Rosseland optical depth scale is indicated by 
the tick marks near the abscissa. All models have identical 
atmospheric parameters (T c g — 2790 K, logg = 5.0, solar 
chemical composition). The dotted line indicates the lo- 
cation of the Schwarzschild boundary in the «mlt = 1.0 
model. See Fig. |l2| for an enlargement of the transition 
region from convective to radiative transport. 
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Fig. 13. Entropy-pressure relation of hydrodynami- 
cal model A-3D (solid) and three standard mixing- 
length model atmospheres (dash-dotted) with oimlt = 
1.0, 1.5, 2.0. The cxmlt = 1-0 model has the largest entropy 
at high pressure while it has the lowest around logP = 5. 
The entropy level of material entering the computational 
box of hydrodynamical model A-3D at the lower boundary 
is indicated by the dashed line. All models have identical 
atmospheric parameters (T e ff = 2790 K, logg = 5.0, solar 
chemical composition). 




4.4 4.6 4.8 5.0 5.2 5.4 5.6 
loflio p B a S [dyn/cnf] 



Fig. 12. Temperature-pressure relations like in Fig. [Tl] but 
focusing on the transition region from convective to ra- 
diative energy transport. In this region the «mlt = 1-0 
model is the coolest, o:mlt = 2.0 the hottest mixing-length 
model. 

transport in LHD. In LHD effects of the turbulent pres- 
sure due to convective velocities computed from MLT are 
neglected. This approximation is well justified since the 
turbulent pressure — as found in the RHD simulations — 
nowhere exceeds 0.5 % of the gas pressure. 

For the hydrodynamical model we have plotted in fact 
six different mean stratifications in the afore mentioned 
figures, each one a horizontal and temporal average over a 
250 s time interval — comparable to the turn-over time of 



the convective cells. The statistical variations are so small 
that the profiles are virtually indistinguishable on the scale 
of the plots demonstrating again that we are dealing with 
atmospheres showing little temporal and horizontal fluc- 
tuations. Fig. [ll] shows that the mean thermal structure of 
the hydrodynamical model corresponds closely to those of 
mixing- length models. The mixing- length models them- 
selves do not show a strong dependence on the mixing- 
length. In the deeper, convective layers small horizontal 
entropy fluctuations are sufficient to carry the convective 
flux, and the radiative energy transport is unimportant. 
Here the stratification has to follow essentially an adiabat 
in the pressure-temperature plane. The same reasoning 
holds for the hydrodynamical as well as MLT models. The 
situation is somewhat different in the radiative layers. By 
construction, the mixing-length models have to approach a 
temperature profile corresponding to radiative equilibrium 
conditions. Since the differences among the mixing- length 
models in the convective region are small the radiative 
equilibrium profiles are almost identical as well. In the 
hydrodynamical model the temperature of the higher at- 
mospheric layers (logP < 5.3) is controlled by a competi- 
tion of radiative heating and adiabatic cooling of material 
which overshoots into regions with formally stable entropy 
gradient < (see Fig. |l3|), and, analogously, the ra- 
diative cooling and adiabatic heating of downflowing ma- 
terial. The relative time scales involved in adiabatic and 
radiative temperature changes determine how the temper- 
ature will adjust between the upper extreme — the radia- 
tive equilibrium temperature — and the lower extreme — 
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the temperature achieved during adiabatic expansion of a 
rising mass element in the higher atmosphere. 

Obviously, in the present case, radiative processes 
dominate and the hydrodynamical model closely follows 
the radiative equilibrium temperature. As we have seen in 
the context of the hydrodynamical models adopting grey 
ATLAS 6 opacities this is not necessarily so. In fact, we 
trace back the short radiative time scales (see Fig. |j) to 
the presence of molecular absorption which provides the 
close coupling of the stratification to the radiative equi- 
librium temperature. We speculate here that the situation 
might change significantly when one goes to metal poor 
objects. Convective velocities are likely not to be reduced 
as dramatically as the atmospheric opacities, shifting the 
relative importance towards adiabatic cooling. This might 
lead to a significant deviation from radiative equilibrium 
in the atmospheric layers of such objects. Indeed, for metal 
poor dwarfs at about solar effective temperatures such ef- 
fects are found in hydrodynamical models (e.g. Aspluiuj 
& Gar cia Perez 200l| ). 

Quantitatively, there are small differences in the ther- 
mal structure between the hydrodynamical model and the 
mixing-length models. As one can expect the differences 
are most pronounced in the transition region between 
convectively and radiatively dominated energy transport 
where the detailed balance between both processes decides 
about the resulting temperature profile. Figures |lj and [l3] 
show that the hydrodynamical model has a mean ther- 
mal structure which is close to the profiles of the mixing- 
length models but cannot be matched exactly. This cer- 
tainly comes not as a big surprise when one considers the 
simplistic nature of mixing-length theory. However, cer- 
tain aspects of the hydrodynamical model can be matched 
by choosing a mixing-length model from the set which are 
parameterized by omlt- This is most easily done consid- 
ering the entropy since it emphasizes model differences 
(see Fig. |l3|). The entropy of inflowing material of the hy- 
drodynamical model which we interpret as the entropy of 
material located in the deep almost perfectly adiabatically 
stratified part of the stellar envelope is matched best by 
a mixing- length model of aMLT(evo) = 1.5. In global stel- 
lar structure models the entropy of the envelope controls 
to some extent the radius of a star, making cvmlt(6vo) 
most relevant for evolutionary models — hence, our nam- 
ing "evo". For further justification and discussion of our 
interpretation see ISteffen (1993| ) and [Ludwig ct al. (1999| ) . 

If one wants to match the entropy jump of the hy- 
drodynamical model «mlt(As) = 2.1 provides the closest 
fit. The entropy jump is the entropy difference between 
the atmospheric entropy minimum and the asymptotic 
entropy of the deeper layers. The entropy jump provides 
a qualitative measure of the overall available convective 
driving, and in the present case also gives a match to 
the temperature gradient of the deeper photosphere with 
log 



model small deviations from equilibrium conditions deter- 
mine the mixing-length parameters. This is a challenging 
situation for any numerical scheme, and we consider it 
possible that the actual calibration errors are dominated 
by related systematic uncertainties. 



5.3.1. Importance of the MLT formulation 

Throughout this paper mixing-length parameters are 
given with respect to the formulation of MLT according 
Mihalas (1978| ) (see |Ludwig et al. (1999| ) for details of the 
implementation). The formulation of Mihalas is widely 
used, in particular in the context of stellar atmosphere 
work. But it is by far not the only MLT "dialect" around. 
Another widely used formulation is the original one by 
Bohm-Vitense (1958| ) which is often encountered in the 
context of stellar evolution work. We would like to em- 
phasize that generally the details in the MLT implemen- 
tation have a direct and significant influence on the cal- 
ibration of the mixing-length parameter. This is particu- 
larly true in the present case where convection penetrates 
far into optically thin regions. Typically it is the treat- 
ment of the radiative energy exchange of the convective 
elements in the optically thin limit that distinguishes the 
various MLT formulations. For enabling inter-comparison 
in other context, we like to mention that the MLT formu- 
lation in PHOENIX is similar but slightly different from 
the Mihalas formulation. However, test calculations have 
shown that numerical differences are negligible in the M- 
dwarf regime. Furthermore, in the ATLAS6 stellar atmo- 
sphere code the Mihalas formulation is implemented. MLT 



formulations as given by Kippenhahn & Weigert (1991) 
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> —1. The internal uncertainties in the deter- 
mination of QfMLT(evo) and «mlt(As) = 2.1 are small 
(« ±0.05) as estimated from the statistical fluctuations 
observed in the models. However, in the present M-dwarf 



and Stix (1989) are identical to the Bohm-Vitense formu- 
lation. 

Figure [li] shows a comparison of the entropy profiles 
of standard mixing-length models of the same «mlt = 1-5 
based on the the Mihalas and Bohm-Vitense formulation. 
We concentrate on the entropy jump since it is sensitive 
to the convective transport in the optically thick and thin 
layers. The asymptotic entropy level of the deeper layers is 
largely controlled by the entropy change over the optically 
thick layers only. Mihalas and Bohm-Vitense formulation 
differ little here, implying the rather similar outcome. The 
formulations differ significantly in the optically thin layers. 
Clearly, the entropy jump comes out to be almost a factor 
of 2 larger in the Bohm-Vitense formulation! To reduce 
the entropy jump to a level comparable as obtained from 
the Mihalas formulation, one had to increase the mixing- 
length parameter noticeably. Again, this emphasizes that 
Q:mlt is only well defined with reference to a specific MLT 
formulation. 

The Bohm-Vitense formulation assumes optically thick 
convective elements throughout. Their radiative energy 
exchange is described in diffusion approximation. It might 
be surprising that the entropy jump obtained from the 
Bohm-Vitense formulation is larger than from the Mihalas 
formulation which distinguishes between optically thick 
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Fig. 14. Entropy-pressure relation of hydrodynamical 
model A-3D (solid) and standard mixing-length model 
atmospheres with omlt = 1-5 based on two different 



formulations of MLT : [Mihalas (1978[ ) (dot-dashed) and 
Bohm-Vitense (1958| ) (dotted). Note, the significantly 
different entropy jumps in the MLT models employing the 
same mixing-length parameter. The entropy level of ma- 
terial entering the computational box of hydrodynamical 
model A-3D at the lower boundary is indicated by the 
dashed line. All models have identical atmospheric pa- 
rameters (Teg = 2790 K, logg = 5.0, solar chemical com- 
position). 



and thin elements. The optically thin elements loose (or 
gain) energy in a non-diffusive fashion. A larger entropy 
jump implies larger radiative losses of the convective ele- 
ments. Intuitively, one would expect that a diffusive trans- 
port is rather inefficient, implying a smaller energ y ex- 
change and smaller entropy jump. However, Eq. ( A. 3 ) 
shows that the situation is reversed: the thermal adjust- 
ment time of convective elements becomes independent of 
their optical thickness in the optically thin limit, while the 
diffusion approximation predicts (inaccurately) a decrease 
of the adjustment time with decreasing optical thickness. 
Since the optical thickness becomes very small in the op- 
tically thin limit, the diffusion approximation underesti- 
mates the thermal relaxation time, and overestimates the 
radiative losses which leads to a larger entropy jump. 

While it appears at first glance arbitrary which MLT 
formulation one chooses, Fig. [l4| also shows that in the 
present case the formulation of Mihalas is perhaps supe- 
rior to the formulation of Bohm-Vitense. The shape of the 
entropy profile of the Mihalas model resembles the result 
of the hydrodynamical model more closely. This comes 
at no surprise since the formulation of Mihalas considers 
optically thick and thin convective elements — i.e. bet- 
ter represents the actual situation. Of course, our choice 
of relating our «mlt to the formulation of Mihalas was 
motivated by this property. 



5.3.2. Comparison to the Sun 

It is interesting to compare equivalent mixing-length pa- 
rameters derived for the thermal structure of the hydro- 
dynamical models with other ones found in hotter objects, 
especially the Sun. Using the MLT formulation of Bohm- 



Vitense, Ludwig et al. (1999) found in a study of the 
convective efficiency based on 2D hydrodynamical models 
c*mlt = 1.6 for the Sun which has likely to be corrected 
to cvjvilt = 1-8 if one considers 3D models. In the Sun no 
distinction needed to be made between fitting the entropy 
jump or the asymptotic entropy. The solar mixing- length 
models match the photospheric entropy minimum found 
in the hydrodynamical models quite well. This is a con- 
sequence of the fact that convection does not reach very 
far up into optically thin regions. Together with the rel- 
atively large entropy jump which makes the exact value 
of the entropy minimum less important this leads to the 
"one «mlt fits both" situation. Contrary to the situation 
in M-dwarfs differences between the Bohm-Vitense's and 
Mihalas' MLT formulation are small in the Sun, again, a 
consequence of the confinement of convection to optically 
thick regions. This allows us to inter-compare the mixing- 
length parameter obtained by Ludwig et al. (1999| ) for the 
Sun with the values derived here. 

The mixing- length parameter of 1.8 for the Sun and 
1.5 (matching the asymptotic entropy) or 2.1 (matching 
the entropy jump) for the M-dwarf are not vastly differ- 
ent, the solar value even falls within the range spanned by 
the M-dwarf values. In MLT the convective flux scales as 
the square of «mlt, hence we are considering a flux vari- 
ation within a factor of 2 only. Convection is efficient in 
M-dwarfs in the sense that the whole envelope is almost 
adiabatically stratified. The convective energy transport 
process itself does not work much differently than in the 
solar case — taking the similar mixing-length parameters 
as a measure. The adiabatic behavior follows mostly from 
the reduction of the transported flux which can be ac- 
commodated with smaller horizontal entropy fluctuations. 
This ultimately leads to the smaller entropy jump which 
amounts to As = 7 x 10 6 erg/ g/K, or only s=s 4 % of the so- 
lar entropy jump. This is remarkable since as evident from 
Fig. |l3| about half of the entropy jump happens to take 
place in the optically thin region which is qualitatively 
different from the Sun. Moreover, the MLT captures sur- 
prisingly well the convective transport properties under 
such conditions. 

We note that an extrapolation of the entropy jump as 
found by |Ludwig et al. (1999| ) (their Fig. 4) for F-, G-, 
and K-dwarfs leads to a significant overestimation of the 
entropy jump for M-dwarfs. This means that the trend in 
the entropy jump found in hotter objects changes at some 
point when moving towards lower T e g. This, again, is plau- 
sible because the dominant atmospheric opacity source as 
well as contributors to the specific heat change in the M- 
dwarf regime. 
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Table 2. Synthetic broad band colors V-I, R-I in the 
Cousins system, as well as J-K, and J-H in the CIT system 
for the models displayed in Fig. [H]. The effective temper- 
ature T e ' ff corresponds to the total flux as computed by 
the detailed spectral synthesis. 



Model 




2rff 


V-I 


R-I 


J-K 


J-H 






[K] 


[mag] 


[mag] 


[mag] 


[mag] 
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= 1.0 


2743 


3.733 


2.010 


0.739 


0.449 


LHD, «MLT = 


= 1.5 


2744 


3.701 


1.996 


0.746 


0.452 


LHD, o-mlt = 


= 2.0 


2745 


3.686 


1.989 


0.749 


0.454 



5.4. Synthetic broad band colors 

Table |^ provides synthetic optical and near-infrared col- 
ors for the models discussed before and shown in Fig. [Tl| . 
The colors of the hydrodynamical model have been cal- 
culated from the average pressure-temperature structure. 
In the calculation of the colors the radiation field has 
been sampled with a resolution of 2 A, i.e. the treat- 
ment of frequency dependence is different and more de- 
tailed than the OBM used in the models proper. Here 
we are concerned with the relative behavior among the 
models. For the V-I and J-K colors appreciable discrepan- 
cies exist between observations and model predictions by 
PHOENIX ("NextGen" model grid) in the r cff -range 
3600 K > T off > 2300 K flChabrier et al. 2000| ; |Allard 
et al. 2 001 ). In V-I the offset amounts to ~ 0.2 — 0.3 mag 
in the sense that the NextGen models are too blue. In 
J-K the offset amounts to ~ 0.2 — 0.5 mag depending on 



the source of water line opacity data used in the models. 
The treatment of the convective energy transport in the 
MLT framework — employing a standard Qmlt of 1.0 — 
was suspected to be one possible reason for the discrep- 
ancy. Table || shows that a detailed treatment of convec- 
tion does not remedy the situation. The difference between 
the V-I color of the hydrodynamical model to the LHD 
o:mlt = 1-0 amounts to 82mmag only. Moreover, the hy- 
drodynamical model is even bluer than the mixing-length 
model, i.e. when assuming the same differential behav- 
ior for models employing detailed radiative transfer, the 
discrepancy between theory and observations gets even 
slightly worse. As expected from the thermal structure 
c^mlt = 2.0 gives the closest match in color among the 
three LHD models. 



Our results strengthen those of Allard et al. (2000) 
and Chabrier et al. (2000) who explain most of the visual 
discrepancy in terms of missing opacities. We would also 
like to emphasize that an optical color mismatch of this 
magnitude does not indicate a major shortcoming of the 
PHOENIX models. The V-band, which is likely to be 
responsible for the color mismatch, is formed high up in 
the atmosphere. The radiative energy flux in the V-band is 
a small fraction of the total radiative output of the star. 



Hence, the problem of the optical colors is related to a 
part of the atmosphere that is of minor importance for 
the overall energetics. 

Similarly, the discrepancies between observed and pre- 
dicted near-infrared colors cannot be traced back to an 
inadequate treatment of convection with MLT. Possible 
differences between MLT and hydrodynamical models 
amount to 10 mmag only. 

Table || allows us to address a side point, namely the 
question of how the hydrodynamical model A-3D and the 
mixing-length models shown in Fig. [ll] can have the same 
effective temperature (2790 K) despite the fact that the 
hydrodynamical model is almost always noticebly hotter 
than the mixing-length models. Indeed, when calculating 
the total flux within the detailed spectral synthesis the 
corresponding effective temperature T c ' ff of the hydrody- 
namical model appears to be slightly (« 10 K) hotter than 
the one of the mixing-length models — with an offset to 
the nominal T e g due to the more detailed treatment of 
the radiative transfer. We interpret this as an effect of the 
small remaining non-linearities when interchanging the av- 
eraging of the temperature-pressure structure with the so- 
lution of the radiative transfer. Strictly speaking, in order 
to determine the colors of the hydrodynamical model one 
must perform the detailed calculation of the 3D radia- 
tion field for many instances in time and had to average 
subsequently horizontally and temporally. However, in the 
present context, the smallness of the effects did not appear 
worth the dramatically increased computational effort. 

5.5. Mean velocity and corresponding effective 
mixing-length parameter 

In the previous sections we discussed mixing-length pa- 
rameters associated with the thermal properties of our 
hydrodynamical models. We now want to turn to a dy- 
namical feature — the velocity profile. Figure [l5] shows 
a comparison of velocities of hydrodynamical model A- 
3D with convective velocities from mixing- length models. 
In the case of the hydrodynamical model RMS fluctua- 
tions of the vertical velocity, for the mixing- length models 
the MLT convective velocity are displayed. Six different 
temporal windows (each 250 s in length) have been used 
for averaging the hydrodynamical data, showing noticable 
variations of the mean velocity particularly in the higher 
atmospheric layers. The velocity plateau in the uppermost 
layers is an artifact of the upper boundary condition which 
stipulates the same velocity at the boundary and the grid 
points next to the boundary. In order to match the convec- 
tive velocities in the convectively unstable layers we find 
that ttMLT(") = 3.5 gives a reasonable match around the 
velocity maximum. The hydrodynamical model has the 
tendency to retain higher velocities towards smaller at- 
mospheric pressure as compared to the MLT models. The 
behavior of the velocities towards high pressure shows a 
faster decline than the MLT models. The deeper hydrody- 
namical model B-3D (not shown) exhibits a similar veloc- 
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Fig. 15. Velocity-pressure relation of hydrodynamical 
model A-3D (solid) and five standard mixing-length 
model atmospheres (dash-dotted) with omlt = 1-0 to 
3.5 in steps of 0.5. The velocities of the mixing-length 
models increase monotonically with «mlt- For the hydro- 
dynamical model RMS vertical velocities from six different 
temporal windows are plotted. A logarithmic Rosseland 
optical depth scale is indicated by the tick marks near the 
abscissa. All models have identical atmospheric parame- 
ters (T c s — 2790 K, logg = 5.0, solar chemical composi- 
tion). 



ity maximum as model A-3D but a decline which is not 
as fast. We conclude that the decline in model A-3D is to 
some extend influenced by the lower boundary condition 
and an asymptotic behavior towards high pressure closer 
to MLT predictions is likely. 

The hydrodynamical model shows an appreciable 
amount of overshoot into the formally stable atmospheric 
regions. Of course, by construction, standard MLT is un- 
able to make any predictions about velocity amplitudes 
there. Two different components contribute to the veloci- 
ties in the higher atmosphere: one contribution are waves 
excited by the stochastical fluid motions in the deeper 
layers and travelling upwards, the other contributions are 
advective motions overshooting into the stably stratified 
layers (see also Ludwig &: Nordlund 2000). 



Qualitatively, we expect that velocity fluctuations v' 
associated with advective motions decay with increas- 
ing distance from the Schwarzschild boundary, roughly 
v' oc Pf where / is a positive parameter of order unity. 
For undamped sound waves we would expect a height de- 
pendence of roughly v' oc P~2. This means that the rel- 
ative importance of both contributions shifts from advec- 
tion dominated motions to wave dominated motions with 
increasing height. If the wave motions are not strongly 
damped we would ultimately expect an increase of the 
velocity fluctuations with height. Indeed, according to 
Fig. [l^ this is what we observe in our model. The wave 
motions mostly arise from standing waves — i.e. excited 
acoustic eigenmodes of the computational box. The de- 
tailed run of the velocity in the overshoot layers reflects 



the location of nodes and anti-nodes of the eigenmodes in 
vertical direction which explains the non-monotonic be- 
havior of the atmospheric velocity amplitude. 

The spectrum and the structure of the eigenmodes are 
to some extent determined by the geometry of the com- 
putational domain. This means that the velocities in the 
higher atmosphere where the velocity amplitude is wave 
dominated are model dependent and should not taken for 
real without skepticism. In the next sections we will even 
argue that the velocity amplitudes of the acoustic modes 
observed in our particular models are probably strongly 
influenced by the numerics, and in this sense are mostly 
artificial. In order to study the transport properties of the 
atmospheric velocity field we will try to filter out as much 
of the acoustic contribution to the velocity field as possi- 
ble. 



5.6. Atmospheric overshoot and transport time scales 

Stellar atmospheres around T c s rj 2800 K are too hot to 
allow for a significant formation of dust grains. However, 
already at slightly cooler effective temperatures, grain for- 
mation sets in and dust grains become major opacity con- 
tributors, i.e. an important factor in determining the ther- 
mal structure of the atmosphere. In fact, the spectral en- 
ergy distribution in the range of effective temperatures 
2500 K < T < 1500 K is crucially linked to the distribution 



of dust grains in the atmosphere (Allard et al. 2000). The 
amount of dust which is present is determined by chemi- 
cal condensation and evaporation processes as sources and 
sinks, as well as macroscopic transport processes which 
carry dust grains away from their sites of formation. In 
M-dwarf atmospheres the transport is dominated by two 
opposing processes: gravitational settling of dust grains 
and their mixing due to the presence of velocity fields, ei- 
ther related to convection or global circulations induced 
by rotation (T. Guillot 2000, private communication). In 
our models, no dust formation takes place but we nev- 
ertheless find it worthwhile to give a characterization of 
the atmospheric mixing due to convection and convective 
overshoot. This can give at least a first order approxima- 
tion of how convective mixing might operate when dust is 
actually present. The basic question which we want to ad- 
dress is whether convective overshoot can provide enough 
mixing to prevent dust grains from being completely re- 
moved from the atmosphere by gravitational settling. 

Formally, we are interested in a statistical representa- 
tion of the mean transport properties of the convective 
velocity field in vertical direction. The horizontal advec- 
tion of dust grains by the convective velocity field probably 
produces horizontal inhomogeneities in the dust distribu- 
tion. We neglect these here since i) we are targeting at 
the application of our results in standard ID stellar at- 
mosphere models, ii) the horizontal inhomogeneities are 
small scale (the size of a convective cell), i.e. hardly ob- 
servable and iii) the uncertainties in our understanding of 
the dust formation process itself perhaps limits the achiev- 
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Fig. 16. Mass exchange frequency as a function of pres- 
sure for various degrees of subsonic filtering. The solid 
curves from top to bottom correspond to no filtering, 
^phasc < 1-0, 0.5, 0.25 km/s, respectively. The dashed line 
is an extrapolation towards the asymptotic behaviour cor- 
responding to a complete removal of acoustic contribu- 
tions. The dashed-dotted line depicts the dependence if 
one assumes a velocity field calculated from linear theory 
A logarithmic Rosseland optical depth scale is indicated 
by the tick marks near the abscissa. The dotted line in- 
dicates the location of the Schwarzschild boundary. 



able accuracy anyway. In view of the last point we present 
a proxy of the convective mixing only, without trying to 
derive a detailed statistical description of the transport 
properties of the convective velocity field, i.e. extracting 



its effective transport coefficients (c.f. Miesch et al. 2000, 
as an example of a more stringent treatment). 

Figure |l6| shows the mass exchange frequency /exchange 
which provides an approximate measure of the convective 
mixing. We define this quantity as 



,/ exchanged ) 



m co \){z) 



(2) 



where F^ ss is the upward directed mass flux, m co \ the 
mass column density, and z the geometrical height. (.} 
denotes temporal and horizontal averages. Note that the 
average total mass flux {F mlLSS } vanishes 



(F mass ) = (pu z ) = (FZ SS + F*Z£) = 



(3) 



p denotes the density, u z the vertical velocity component, 
/exchange is the frequency with which the mass above a 
certain height z is replaced by the flux of fresh material 
from below. Mass conservation and quasi-steadiness de- 
mand that the upward directed mass flux from below is 
compensated by a downward directed flux from above. As 
long as the flow field is sufficiently disordered this mass 
exchange should also go hand in hand with mixing of ma- 
terial; in this way, /exchange can be interpreted as mixing 
frequency of atmospheric layers above a certain reference 
height z. The assumption of a disordered velocity field 
was checked qualitatively by studying the behaviour of 



passive tracer particles being advected with the flow. The 
definition of /exchange was motivated from our expectation 
that the amplitude of the velocity field would strongly de- 
cline with distance from the convectively unstable zone. 



Frcytag et al. (1996) find in hydrodynamical models that 



the velocity amplitude due to convective overshoot de- 
clines exponentially with distance from the Schwarzschild 
boundary, and promote this result as a generic feature of 
overshoot. In such a situation, the mass flux in or out of 
a control volume is dominated by the flux through the 
surface closest to the convectively unstable region. This is 
the reason why we consider just one surface (or height) 
through which mass is transported. 

As mentioned previously, the atmospheric velocity field 
has strong contributions of acoustic waves generated in the 
deeper layers of the models. From the comparison (not 
shown) of the hydrodynamical models A-3D and B-3D — 
which has a slightly greater extent in depth — we con- 
cluded that instabilities associated with the lower bound- 
ary condition contribute to the excitation of waves, i.e. 
the waves are mostly a numerical artifact. We could not 
trace down the responsable process, but the behavior was 
apparently related to the high degree of adiabaticity the 
stratification exhibits close to the lower boundary. In the 
Sun the boundary condition operated stably. Besides the 
numerical findings, from a purely theoretical point of view 
one would expect that the generation of acoustic energy 
is rather inefficient in a flow with as low a Mach num- 
ber (few percent) as in the present case. Finally, it turns 
out that the Richardson number associated with the hor- 
izontal shear introduced by the acoustic modes in the at- 
mosphere is not low enough to render the shear flow un- 
stable. No small scale turbulence is induced by the waves 
and their oscillatory velocity fields do not contribute sig- 
nificantly to the mixing of material. Hence, we concluded 
that the wave contribution to the velocity field should be 
discarded if one wants to measure its mixing efficiency. 

We therefore removed the contribution of waves before 
evaluating the average mass flux by subsonic filtering — a 
technique invented in the context of solar observations for 
cleaning images from "noise" stemming from the 5 minute 
oscillations ( Title et al. 198S ). In short, one considers a 
time sequence of images and removes features with hori- 
zontal phase speeds greater than a prescribed threshold. 
This is achieved by Fourier filtering of the spatial-temporal 
data in the k-ui domain. In our case we worked with 150 
snapshots of the flow field sampled at 10 s time intervals. 
The total length of the sequence was long in comparison 
to the periods of the excited acoustic modes. For every 
depth layer we performed a 3D Fourier analysis (one tem- 
poral, two spatial dimensions) of the vertical mass flux. 
We selected velocity thresholds 



Wphasc = jT7 = 1-0, 0.5, 0.25 km/s. 



(4) 



With these specific choices of v p hase we tried to remove 
acoustic features which had phase speeds in the order of 
the sound speed without affecting the slow convective mo- 
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tions which have maximum velocities of about 0.25 km/s 
(see Fig. |l5|). Note that the technique accounts only for 
the horizontal component of the phase speed. This means 
that e.g. purely radial oscillations could not be removed 
by this procedure. An improvement could, in principle, ex- 
tend the filtering technique to four dimensions (one tem- 
poral, three spatial dimensions). But this would require 
a much higher computational effort, in particular since 
one would first have to determine a suitable set of basis 
functions for the representation of features in the vertical 
direction, where simple harmonic functions do not form 
a proper set of eigenfunctions. Luckily, the standard sub- 
sonic filtering worked well enough in our case so that we 
did not have to pursue this issue further. 

Figure [l6] shows the resulting mass exchange fre- 
quency /exchange as a function of height (expressed in 
terms of the mean pressure) for our hydrodynamical 
model A-3D. Lowering the threshold v p hasc leads to a sig- 
nificant decrease of the mass exchange frequency in the at- 
mospheric layers while the deeper layers remain essentially 
unaffected by the subsonic filtering. This again shows that 
the velocity field of the higher atmospheric layers stems 
primarily from wave motions. Only the lowest threshold of 
0.25 km/s starts to affect the deeper layers since also con- 
vective flow features reach such velocities. The dashed line 
in Fig. [l6| indicates our extrapolation of what to expect in 
case of a perfect removal of the acoustic contribution. Its 
slope is given by the envelope of the sequence of curves for 
increasing degree of subsonic filtering. We want to empha- 
size that by the very nature of the problem the filtering 
procedure does not apply a small correction to the sig- 
nal but removes even most of the original signal in some 
layers. Wave motions almost inevitably dominate the flow 
velocity above a certain atmospheric level. Large correc- 
tions as exhibited in Fig. |l^ are to be expected, but do not 
pose particular problems so that we are confident about 
the filtering procedure. 

As an additional check of our approach we calculated 
for an LHD mixing- length model with bmlt = 2.0 linear 
eigenmodes describing the growth of a perturbation un- 
der the action of the convective instability, omlt = 2 was 
chosen since the model gave among the available mixing- 
length models the closest match to the thermal struc- 
ture of the hydrodynamical model. As horizontal wave- 
length of the pertubation we chose the horizontal extent 
of model A-3D of 250 km, i.e. the largest horizontal wave- 
length which could be accomodated in this hydrodynam- 
ical model. Freytag et al. (1996| ) showed that the velocity 
field of convective eigenmodes gives a reasonable repre- 
sentation of the velocity field in the overshoot region as 
long as the horizontal wavelength is chosen in the vicinity 
of the actually present convective scales. The dash-dotted 
line in Fig. [l6| shows the shape of /exchange associated with 
the linear eigenmode. The velocity amplitude was scaled 
to match the maximum of /exchange in the deeper layers of 
the hydrodynamical model. Besides a systematic shift the 
functional behavior predicted by the mode is quite similar 
to the one of the hydrodynamical simulation. Moreover, 



the time scale for the growth of the perturbation comes 
out to 120 s which is of the same order as the convec- 
tive turn-over time scale. The similar rate of decline of 
/exchange in the atmospheric layers strengthens our con- 
fidence in the robustness of our extrapolation procedure. 
However, changing the horizontal wavelength of the per- 
turbation within reasonable limits as well as experiment- 
ing with different background models did not improve the 
correspondence between linear modes and non-linear hy- 
drodynamical model beyond the quality shown in Fig. 

Strictly speaking, the mixing frequency displayed in 
Fig. |l^ is a lower limit since we removed all contributions 
associated with waves. However, the arguments given be- 
fore make it very likely that the contributions of wave 
motions to the mixing are negligible: the Mach number of 
the flow is low, and the weak waves motions which may be 
excited are associated with little turbulence. Therefore we 
believe that the derived mixing efficiency is close to what 
is encountered in a real M-dwarf atmosphere. 

What are the possible consequences of the exponential 
decline of the mixing time scale due to convective over- 
shoot? Following the work of Rossow (1978), one arrives 
at an estimate of about 10 4 s for the typical time scale of 
dust sedimentation in M-dwarf atmospheres. According to 
our hydrodynamical models the convective mixing could 
counteract the sedimentation up to a layer of about two 
pressure scale heights above the formally convectively un- 
stable layers. The convective turn-over time scale is of the 
order of 10 2 s. This is slow enough to allow for dust nucle- 
ation in the upper regions of the convective envelopes of 
M-dwarfs — provided that the thermodynamic conditions 
allow for nucleation in the first place. From these consider- 
ations we would expect that dust clouds in M-dwarf (and 
brown dwarf) atmospheres are confined to layers in the 
vicinity of the upper Schwarzschild boundary of the con- 
vective envel 



6. Concluding remarks and outlook 

We used elaborate 2D and 3D radiation-hydrodynamics 
simulations to study properties of convection on the sur- 
face of a prototypical late M-dwarf (T e g « 2800 K, log p = 
5.0, solar chemical composition). Despite the significant 
differences in the physical conditions encountered in the 
solar and an M-dwarf atmosphere we obtained the strik- 
ing result that M-dwarf granulation does not look qual- 
itatively different from what is familiar from the Sun 
(see Fig. Quantitative differences (intensity contrast 
1.1%, horizontal scales « 80 km, maximum RMS veloci- 
ties « 240 m/s, convective turn-over time scale ~ 100 s) re- 
main within the expectations derived from mixing-length 
theory. 

Connected to this basic finding is the — for practical 
purposes — perhaps most important result that the tem- 
perature structure of the higher atmospheric layers is de- 
termined by the condition of radiative equilibrium, and is 
not very much affected by processes usually not accounted 
for in standard stellar atmosphere models. Convective 
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overshoot as well as energy transport by waves do not sig- 
nificantly affect the temperature structure outside of the 
region of convective instability. We expect that this find- 
ing also holds for main-sequence objects of higher effective 
temperature where radiation becomes relatively more im- 
portant. 

Answering the questions posed in the introduction we 
confirm that in late M-dwarfs mixing- length theory can be 
applied to obtain a realistic description of the convective 
energy transport in a ID stellar atmosphere code, pro- 
vided one can remove uncertainties related to the choice 
of the mixing-length parameter. If an independent cali- 
bration is available we expect that ID atmosphere models 
allow reasonable accurate predictions (on a level as dis- 
played in Fig. 13) of the atmospheric temperature struc- 
ture and ultimately the stellar spectrum to be made. 
Effects on the stellar spectrum related to horizontal tem- 
perature inhomogeneities are expected to be very small 
in the M-dwarf studied here, primarily due to the small 
horizontal fluctuations of the thermodynamic variables. 
Even in the Sun — with the much higher temperature 
contrast present at its surface — effects related to hori- 
zontal temperature inhomogeneities are rather subtle (cf. 
Stcffcn fc Ludwig 1999| ). We expect that ID stellar atmo- 



sphere models provide an acceptable overall approxima- 
tion to the spectrum of main-sequence objects between 
the Sun and late M-dwarfs. Only if one demands for a 
precision exceeding commonly adopted levels or desires to 
study effects in principle not included in standard model 
atmospheres (e.g. spectral line shifts and asymmetries), 
one has to go to more sophisticated modeling. We em- 
phasize that this refers to objects at solar metallicity. For 
metal poor objects the situation is clearly different (see 
Asplund ct al. 1999| ). 

A downside of our findings is that we cannot trace 
back the remaining differences between theoretical and ob- 
served colors for M-dwarfs to shortcomings in the treat- 
ment of the convective energy transport. The resolution 
of the problem has to be found elsewhere, as pointed out 
earlier, deficiencies in the molecular opacities are still a 
possible option. 

We are left with the problem of finding an adequate 
value of the mixing- length parameter to obtain a descrip- 
tion of the vertical temperature run in the superadiabatic 
layers. In the present case we find a mixing-length param- 
eter «mlt(As)=2.1 which gives a match to the entropy 
jump and the temperature gradient of the deeper atmo- 
spheric layers (see Fig. |l^). Despite perhaps the best value 
to be employed in stellar atmosphere calculations, for 
global stellar structure models a value of q;mlt(As)=1.5 
would be more appropriate since it ensures to find the cor- 
rect asymptotic entropy in the convective envelope which 
is important for obtaining the correct stellar radius. To 
complete the "zoo" of mixing-length parameters, we get a 
value of q;mlt(i0=3.5 when matching the convective ve- 
locities predicted in our hydrodynamical models. We re- 
iterate that all value s are given with reference to the for- 
mulation of MLT by iMihalas (1978| ). 



The various values of the mixing-length parameter 
point towards the deeper rooted problem that MLT can 
give a reasonable but not exact description of the average 
convective properties. Even calibrating one aspect does 
not ensure the overall correct functional form of, say, the 
temperature profile. We have seen that different formu- 
lations of MLT can give quite different functional depen- 
dencies. They offer the possibility to improve fits beyond 
the quality limited by fitting the mixing-length parame- 
ter only. For late M-dwarfs all this does not matter much 
since the differences of the atmospheric structure for var- 
ious values of «mlt are small. However, we stress that 
this statement refers to cooler M-dwarfs on or close to 
the main-sequence. For pre-main-seq uence (PMS) object s 
the situation is markedly different ( Baraffc ct al. 2002 ). 
There the specific choice of «mlt has a large impact on 
the resulting atmospheric structure. Work is underway to 
extend the present study into this regime which may also 
allow us to find the most suitable MLT formulation. 

A result beyond the scope of classical model atmo- 
spheres is the derivation of a proxy of atmospheric mixing- 
time scales (see Fig [l6]) due to convective overshoot. We 
find an exponential "leaking" of the convective velocity 
field into the formally stably stratified layers. Depending 
on the exact criterion, overshooting extends the efficiently 
mixed regions about 2 pressure scale heights beyond 
the Schwarzschild boundary. We suggest that the mix- 
ing found in the T e g « 2800 K model studied here, takes 
place in an analogous fashion in brown dwarfs, and pro- 
vides the mixing which counteracts dust sedimentation. 
Hydrodynamical models can be used to address this prob- 
lem more directly by performing simulations including the 
formation and transport of dust which we consider as an 
interesting challenge for the future. 

Last but not least we would like to point out the two 
weakest points of our investigation. While the precision 
of the radiative transfer in the hydrodynamical calcula- 
tions is sufficient to address the questions discussed here 
there is certainly room for improvement of the OBM to 
get an even closer agreement with detailed spectral synthe- 
sis calculations. Secondly, the models presented here are 
rather shallow. Improvements in the formulation of the 
lower boundary condition would perhaps allow the use of 
deeper computational domains, and would reduce the in- 
fluence of the specific formulation of the lower boundary 
conditions. 
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Appendix A: Computation of the characteristic 
time scales 

The radiative time scales shown in Fi g. p| were evalu ated 
in Eddington approximation (see e.g. Edwards 1990| ) ac- 
cording 



^rad 



16a X T 3 



3xV 



(A.l) 



(c p denotes the specific heat at constant pressure, a 
Stcfan-Boltzmann's constant, x opacity, T temperature, 
p mass density, k wavenumber of the disturbance) . We as- 
sumed that the geometrical size of a thermal disturbance 
is the local pressure scale height Hp, hence 
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Introducing the optical size of the disturbance Tdi S = 



XpHp Eq. (A.l) can be written as 



^rad — 



1 



4tt 2 



(A.3) 



Relation (A.3) shows that in the optically thin limit Tdi S <C 
1 the radiative relaxation time is independent of Tdi S , and 
in the optically thick limit Tdi S S> 1 scales quadratically 



with Tdis- Equations (A.l) and (A.3) describe the decay 
of a Fourier mode in an infinite isothermal atmosphere, 
i.e. make very idealizing assumptions. Our choice of of the 
local pressure scale height as size of the perturbation is 
somewhat arbitrary. Hence, the reader should take i ra d as 
depicted in Fig. |] as order of magnitude estimates. 

The (adiabatic) Brunt- Vaisala period £bv was calcu- 
lated according 



2n 



BV 



: , with w|v = S- (V ad - V) . (A.4) 
tip 



(5 = — ( j^r j denotes the thermal expansion coefficient 

at constant pressure, g gravitational acceleration, V ad the 
adiabatic, and V the actual temperature gradient of the 
atmosphere). The local Kelvin-Helmholtz time scale £kh 
was evaluated according 



Pc p T 



(A.5) 



(P denotes the gas pressure) under the assumption that 
the overall energy content of a mass column is dominated 
by its thermal energy content. 
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